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We present a maximum entropy thermodynamics to analyze the state space of a small 
system interacting with a large bath. In small systems, the fluctuations around the mean 
values of observables are not negligible and the probability P{i) of the state space {i} of the 
system cannot be described by a unique set of parameters C that characterize the interaction 
of the system with the surrounding bath. We employ a superstatistical approach: The 
probability distribution P{i) for the phase space {i} is expressed as a marginal distribution 
summed over the variation in The joint distribution P{i, ^) is estimated by maximizing 
its entropy. 

We test the development on a simple harmonic oscillator strongly coupled to a bath of 
Lennard- Jones particles. The estimated distribution P{r) of the position r of the oscillator 
depends on the information that is used to construct it and not all measurements have equiv- 
alent predictive power. Moreover, the traditional canonical ensemble distribution emerges as 
a limiting case of a much richer class of niaxEnt distributions. Future directions and other 
connections with traditional statistical mechanics are discussed. 
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I. 



INTRODUCTION 



Thermodynamics is a science of large systems. In thermodynamic systems, the fluctuations 
around the mean value of observables are negligible compared to the mean and the statistical 
mechanical estimation of the probability distribution P{i) of the states {i} of the system depends 
only on the global functions of the states such as total energy E, total volume V, total magnetization 
M, etc (jHEj)- The maximum entropy (maxEnt) interpretation views statistical mechanics as an 
inference problem — the distribution P{i) of states {i} is estimated from the limited available 
knowledge of the system. Briefly, the maxEnt program of estimating probabilities P{i) of states 
{i} of a system coupled to a bath involves maximizing the entropy function S[P{i)] subject to 
constraining the values of certain experimentally known observables of the system. For example, if 
{Xk} are the mean values of the fluctuating observables {X^}, then the probabilities P{i) of states 
{i} are estimated by maximizing the function in Eq. [l] ([H El i]) . 



{C^} and 7 are Lagrange multipliers that ensure that the imposed constraints are satisfied and that 
the probabilities are normalized respectively. The entropy is a non-negative convex function of the 
probabilities and is usually defined as ([51 IH) 



The estimated distribution P{i\(^) is parametrized by a unique set of Lagrange multipliers that 
ensure that the ensembles averages 




(1) 




(2) 




(3) 



are equal to the experimentally observed ones (i.e. {{Xk) = Xk}) and is given by 





Above 




(5) 
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is the partition function. The Lagrange multipHers are determined by solving 



— ^7t — " 



Owing to the law of large numbers, for thermodynamically large systems, the relative spread 
around the mean value is negligible compared to the mean value for any obserable quantity (jll 
[3 |8]). Traditional statistical mechanics enjoys great success in predicting the behavior of such 
systems. This may occlude its inference aspect: for thermodynamic systems, the mean values of 
the observables {X^} are sufficient to correctly predict the system behavior. Additional information 
about their higher moments does not result in an improved predictive statistical mechanics (4). In 
fact, the higher moments {XJ^) (n > 1) can be correctly estimated from the partition function ([21 [7|). 

The key feature of real systems that make statistical mechanics a success is the sharply peaked 
distribution estimated by Eq. [T]([T]). Sufficient (but not necessary) conditions for the same include 
a) locality of interactions, b) extensivity, and c) large size of systems. The accuracy of traditional 
statistical mechanics may not be carried over to systems that violate either of these requirements. 
Modifications to the maxEnt program are necessary to describe such systems. Some examples of 
such modifications include making the entropy functional non-extensive ((Qj) or including higher mo- 
ments of the observables as constraints (jlOp . Small systems do not satisfy these requirements (jTl IH)) . 
Thus, according to the maxEnt interpretation, unlike for a thermodynamic system, the predictions 
about a small system should depend on the information that is used to estimate the distribution 
P{i). 

In this article, we present a maxEnt generalization of traditional statistical mechanics towards 
predictive statistical mechanics that is applicable to both small systems and large systems. We 
elucidate the approach with an example of a harmonic oscillator coupled to a bath of Lennard- 
Jones particles. We choose the harmonic oscillator in this initial study since it is one of the few 
systems whose partition function can be computed analytically. This allows us to clearly illustrate 
the inference aspects of our development. We hope that our method is of general importance to 
the study of thermodynamics of small systems. 

The article is organized as follows. In sec. |llj we describe the theoretical development. In 
sec. 



Ill we analytically work out the and compare some numerical results. In sec. IV we discuss 
connections of our method to traditional statistical mechanics along with possible limitations and 
generalizations. 
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II. THEORY 
A. Maximum entropy thermodynamics 

Imagine a constant volume system coupled to a large bath. The interactions within the system 
and the interactions of the system with the surrounding bath determine the distribution P{i) of 
its states {i}. If the interactions between the system and the bath are weak compared to a) the 
interactions within the system and b) the interactions within the bath, they can be characterized 
by a unique set of parameters where the distribution P{i\^^) is parametrized by . For example, 
the chemical potential /x and temperature T dictate how a fluid within a given macroscopic volume 
V exchanges molecules and energy with its surrounding. 

If the system under consideration is small, the system-bath interactions are non-negligible com- 
pared to the interactions within the system and cannot be characterized by a unique set C"^. For 
example, if the volume V of a ^VT system is comparable to molecular sizes, the so called surface 
effects become important and there no longer exists a unique chemical potential /i that governs the 
average number of particles within the confines of the volume V . In such cases, one must allow 
C, to vary and consequently, the entropy of the joint distribution C) instead of P{i) should be 
maximized ()lip . Thus, the optimization problem involves maximizing 

5[P(i,C)] = -5^^(^'C)iogni,C) 

with the constraints 

i 

= Y,PiC){Xk)c = X-k (7) 

c 

for A; = 1,2, ... and 

{Ym)=^P{i,C)Yrr,{0=Y,n (8) 

for m = 1,2, ...,M. Here, we introduce measurements {Ym} that constrain and dictate the 
variation in the parameters themselves. Note that such measurements are redundant for a ther- 
modynamic system since a thermodynamic system is determined by a unique set of parameters. 
Also note that while {Xk{i)} depend solely on the state space {i}, {Ym{C)} depend on C- 
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To recast the above problem, let us write 



where 



S[P{^,0] = S[PiO] + S[Pm] 



s[pm] = E - E ^(^10 log ^(^10 



(9) 



(10) 



is the conditional entropy of the state space In Eq. 10 we have replaced the summation by 



S'(C) for brevity. S{C^ is the entropy of the system z/ it were to be described by a unique set of 
Lagrange multipliers C. Thus, the objective function that needs to be maximized (including the 
constraints) is (see Eq. [9] and Eq. 10 ) 



+ E"fc 

k 

+ Am 

m 

Summing over \i\ degrees of freedom. 



E^'(i,c)i;n(c) 



^ vr, 



(11) 



^[^'(c)] + E ^(c)s(c) + 7 E - 1 



+ E 



E^(c)(^, 



+ E^^ 



fc/C 



E^(^")^-(o 



Y 



(12) 



Carrying out the maximization 

P{0 



]—-—-exp f 5'(C) - E «fc(^fc)c ~ E ^™^™(C) 



(13) 



-2^({afc}, 

Here, Z{{ak},{\n}) is the partition function. Finally, the marginal probability distribution P{i) 
is given by. 



P{^) 



Z{{ak} 



^— ^E^(^lC)-e^^^('5(C)-E"^'™C-E^'»^-(C) ) • (14) 

^^^^^^1 ? \ k ra j 
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Eq. 14 is the probability distribution P{i) of states {i} of the smaU system which incorporates ah 
of our knowledge about the system viz. the observations {X^} and the deviation of the system 



from the thermodynamic limit i.e. the variability in (see Eq. 12). 



B. Choice of constraints and the thermodynamic limit 
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It is instructive to see if there exists a limiting case of Eq. 
particular value of parameters . This limiting case is the thermodynamic limit of Eq. 
so, let's write = 6uJk a-nd Am = ^Km- Then, 



which reduces to P{i\C^) for a 

To do 



14 



P(i) 



1 



Z{{ak}, {Xm}) 



Y,pm-expls{c)-s 



C 



Wfc {Xk)^ + ^ KmYm (C) 



(15) 



Observe that as (5 — t- oo, only the value(s) of that correspond to local maxima of the 6 dependent 
part of the exponential ujk{Xk)^ + Ylim i^mXmiC) contribute to the summation in Eq. 
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Con- 



sequently, Eq. [14] has a unique thermodynamic limit if above summation has a unique maximum. 
It appears that apart from this requirement, there exists no additional restriction on the choice of 



experimental constraints {1^(0}. Moreover, a unique maximum ensures that Eq. 14 is indeed a 
generalization of Eq. [4] as the latter turns out to be a special case of the former. 



III. HARMONIC OSCILLATOR STRONGLY COUPLED TO A BATH 
A. How does strong coupling manifests itself? 

We will carry out the above program and derive Eq. [14] for a harmonic oscillator coupled to a 
large bath to illustrate, with a concrete example, the effect of strong system-bath coupling on the 
estimated distribution P{i) of states of the system and its subjectivity. 

The internal states {%} for the oscillator are the continuous variable r denoting the deflection 
of the oscillator from its reference. The potential energy U{r) of the oscillator when the spring 
constant is ko is given by, 

U{r) = ^kor'^. 

If the oscillator is weakly coupled to a thermodynamic bath at an inverse temperature f3, statistical 
mechanics estimates the probability distribution P{r\ko,f3) (excluding the volume element Airr^) 
as, 

P{r\ko,l3)(xexp(-^l5korA (16) 
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What happens if the oscihator is strongly coupled to the bath? If the subscript b denotes 
the bath then recognizing that the large system comprising of the oscillator and the bath is ther- 
modynamic in nature, the distribution of states P{r, rt) of the oscillator and the bath is given 

by, 

P{r, n) oc exp ^-/3^A;o?-^ - PU\,{r, rb) - /3C/s(r, rb) 

where C/b is the interaction energy within the bath, Us is the interaction energy of the bath with the 
system (the oscillator), and /? is the inverse temperature. Integrating over rb degrees of freedom, 

P(r) oc (-/3^A:or2-/3</,(r;/3)^ . (17) 

/ 1 2 
p< exp ( —-(3kor 

Thus, if the modulation of the oscillator state space distribution due to Us{r,r];)) cannot be 



neglected (jT2HT5]) . P{r) is no longer estimated by Eq. 16 Moreover, the temperature dependent 
molecular field (p{r; f3) is a non-trivial function of the details of the molecular coupling between 
the bath and the small system and, in general, cannot be analytically determined. The maximum 
entropy development presented above estimates this field from limited experimental information 
about the small system in a maximally noncommittal fashion. Below, we will re-derive the distri- 
bution P(r) and show that the weak coupling picture of traditional statistical mechanics emerges 
as a limiting case of a much richer distribution which is in fact subjective. 

From here ownards, without loss of generality, assume ko = 1. This is equivalent to absorbing 



ko in /3. As above, if we know the mean energy U (r) of the harmonic oscillator coupled weakly 
to a bath, the distribution of the position r is estimated to be P(r|/co,/3), given by Eq. 16, at a 
particular value f3. 



4^2^3/2 



P{r\P) = ^ ^ (18) 



vr 



Here, /3 is the parameter (= (^) of the thermodynamic system and is no longer unique if the 



oscillator-bath coupling is strong (see Eq. 13). Note that /3 also has in it the dimensionless spring 
constant /cq- Thus, a variation in /3 can be interpreted as a variation in /cq and not in the temperature 
T. Mathematically, this treatment is similar to the superstatistical generalization of statistical 
mechanics (1161). 
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B. Case 1: Constraining U{r) and the average entropy S 

The entropy of the oscihator alone is not maximized and a natural measure of the variability 
in /3 is the measured entropy S (17). We introduce the contraint that the ensemble average (see 



Eq. 10) 



{sm = j2msm. (19) 

c 

is equal to the measured entropy S. 

The entropy 5(/3), of the distribution P{r\/3) at a paricular value of /3 is given by 

5(/3)~-^log/3 

upto an additive constant. 



Following Eq. 13 our best estimate of the probability P{f3) of /3 is given by, 



and the marginal distribution 



Pi(/3) ocexpQlog/3-C^ 

(2)1-1 g-i^-A/2^|-l 



Pi(r) = / Pi{l3) ■ P{r\/3)df3 



is given by, 



11 A A+1 A~l A+1 , ^ „^ 

V^r (I - 1) 

Here, K^{x) is the modified Bessel function of the second kind with parameter 7. To understand 



^1 (0 = ^^.A J (21) 



Eq. [20] and Eq. 21 physically, let's write ( = kX and calculate the moments of Eq. ^Tl The first 



two moments are given by. 



As A — > 00, 
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FIG. 1. Left: As A c», the distribution of (3 approaches a diract delta distribution S{f3 — 3k). 

Notice that at small values of A, the Pi{/3) distribution is very broad. Right: Similar to the left panel, as 
A — !■ CX3 the distribution Pi{r) tends to the distribution in Eq. 18 of a harmonic oscillator coupled to a weak 




bath at inverse temperature P ~ 3k (blue curve). Here k ^ (/X. Physically, A measures the strength of the 
coupling between the bath and the small system. 



It is easy to see that as A — ?■ oo, Eq. [21] approaches the canonical ensemble distribution Eq. [18] 



at /3 = 3k. Eq. 18 represents the distribution of the harmonic oscillator only when the coupling 
between the oscillator and the surroudings is weak. Thus the Lagrange multiplier A measures the 
strength of the coupling between small system and the bath while k = ^/ A is the effective inverse 
temperature (or the effective spring constant) of the harmonic oscillator. In Fig. [T] we illustrate 
this graphically. As A — )• oo, P{f3), the distribution of /3 approaches a diract delta distribution 
6{f3 — 3k) implying that the system is described by a single inverse temperature f3. Consequently, 



the distribution P{r\X, ^ = kA) approaches the distribution in Eq. 18 with /3 = 3k. 



C. Constraining U{r) and /3 
Instead constraining S, we can introduce different information about the individual parameters 



in the constrained optimization problem in Eq. 12 For the harmonic oscillator, let us examine 



the consequences of constraining the mean value /3 of the effective inverse temperature /3. The 
derivation is straightforward and we only show the distributions P{13) and P{r), 



m) = (24) 



and the distribution P^ij^ is given by, 



= ■ ^''^ 
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FIG. 2. Left: As C ^ oo, P2{P), the distribution of /3 approaches a diract delta distribution 5{j3 — ^ )■ 
Notice that similar to Eq. [20]at small values of C, the P2{(3) distribution is very broad. Right: Similar to 
the left panel, as C — > oo the distribution P2ir) tends to the distribution in Eq. Ilslof a harmonic oscillator 



coupled to a weak bath at inverse temperature /? — ^ (blue curve). Here k — ^/^. Similar to A in Eq. 
here, C measures the strength of the coupling between the bath and the small system. 
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The moments (r) and (r^) are given by, 




IT 



(26) 



Again, we put = nC aiid take limit C — >■ oo to get, 

2 4/2 4/ 



IT 



(27) 



Thus, similar to Eq. [2TJ Eq. [25] also reduces to the canonical ensemble distribution in the limiting 
case Q ^ oo (see Fig. [2]). The inverse temperature /3 of the oscillator in the thermodynamic limit 
is given by 



(28) 



Briefly, the maxEnt program estimates the probability distribution P{r) (Eq. 21 and Eq. 25) for 



the position r of the harmonic oscillator from the available information about the observable mo- 
ments of the position and information about the system-bath coupling. The maxEnt distributions 
reduce to the traditional statistical mechanical estimate of P{r) in the case where the coupling 
between the oscillator and the bath becomes weak. 
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r r 

FIG. 3. The empirically observed (blaek eircles) and fitted (red and blue) distribution P{r) for the position 
r of the the harmonic oscillator for two different strengths of coupling between the oscillator and the 
surrounding bath (right is somewhat stronger than left, see supplementary materials for details). The 
red curve is the best fit for Eq. [21] and the blue curve is the best fit for Eq. [25) The black curves represent 
the best fit of the canonical ensemble distribution, Eq. [TSj 



D. Numerical simulations 



We test our analytical models, Eq. [21] and Eq. 25 on molecular dynamics simulation of a har- 
monic oscillator at the origin coupled to a gas of Lennard-Jones particles. We choose the coupling 



between the bath and the system such that U{r) « Us{r,r]^) (see Eq. 17) (see supplementary 
materials for details). 

Due to the strong coupling, we expect that the distribution P{r) will be different from P(r|/3). 
Fig. [3] shows that indeed the empirically observed distribution P{r) (black circles) for the two 
simulation systems has an extended tail that cannot be captured by the canonical ensemble distri- 
bution (black line). The two distributions Pi{r) and P2{r) (blue and red respectively) capture the 



empirically observed distribution sufficiently well only from the knowledge of average energy U{r) 
and an additional piece information about the fluctuation in the parameter /3 (see supplementary 



materials for fitting procedure). Moreover, the distribution estimated from the knowledge of U{r) 



and S (Eq. 21 ) is a slightly better description than the distribution estimated from the knowl- 



edge of U{r) and /3 (Eq. 25). This suggests that, at least in the case of the harmonic oscillator, 
the measued entropy is a better quantifier of the deviation of the oscillator system from being 
thermodynamic in nature (I17p than its mean inverse temperature /3 with respect to predicting the 
behavior of the oscillator. 

These results suggest that the distribution of states of a small system can be predicted accurately 
from a few measurements {Xk} and {Ym} of the state space and of the parameters respectively 
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that constrain their fluctuations. Even though it is apparent that the predictions will depend on 



the choice of observables {X^} and especially {Ym} (see Eq. 13 and Eq. 14), it is not very clear as 
to which choices lead to a better prediction. We believe it to be a difficult question and leave it 
for futher investigations (see below). 



IV. CONCLUDING DISCUSSION 



The maximum entropy principle estimates the distribution of states of a system coupled to 
a bath when the system is very large, the bath is very large compared to the system, and the 
interactions between the system and the bath are negligible. In this case, the distribution P{i) 
of states {«} of the system is characterized by a unique set of parameters that characterize the 
bath. If the system under consideration is small, the bath is not characterized by a unique set C 
and instead one must entertain the entire distribution P{C)- 

In this work, we have presented a maximum entropy method to predict the behavior of a 
small system in contact with a large thermodynamic bath with minimal knowledge about the 
small system. Such problems are becoming numerous especially as new technology allows precise 
measurements at small length scales. In the current work, we developed the thermodynamics of a 
small system which depended on the knowledge of a) the mean values of some observables {Xf.} 
and additionally the variability in the parameters C that describe system-bath coupling. 

Here, we have not discussed the possible experimental methods to a) estimate the fluctuations in 
the parameters and b) measure the variables {^^(C)} that constrain the variation in the parameters. 
Fortunately for the harmonic oscillator, for operational purposes, the specifically chosen variables 
{Ym{C)} could be mapped to suitable observables of the state space {z},i.e. we could find Xmii) 
such that (see supplementary materials), 

Y„^{C) = Y,P{i\OXm{i). (29) 

i 

We leave it for further studies to investigate the generality of the above mapping. 

The current work also shows that the maxEnt framework developed in (17) for non-equilibrium 
systems is also applicable to small systems. Moreover, we suspect that the key findings in the 
current work viz. the subjective nature of estimated probability distributions for small systems and 
their objective traditional statistical mechanical limit both are features of that framework. 

One criticism of the maximum entropy interpretation of statitical mechanics is that it does not 
lead to predictions that are otherwise inaccessible to traditional statistical mechanical methods. 
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To the best of our knowledge this is the first work that clearly highlights the inference aspects 
of predictive statistical mechanics ([1]) and makes predictions that are difficult to come by via 
standard statistical mechanical techniques unless one knows the details of the interaction of the 
small system with the bath. The harmonic oscillator allows us to work analytically and we show 
that the traditional canonical ensemble is a limiting case of a much richer probability distribution 
that is in fact subjective in nature: it depends on our knowledge of the system (see the tail region of 
the distributions in Fig. [s]). The subjectivity of the estimated distribution is sometimes considered 
to be a weakness of the foundations of maxEnt interpretation of statistical mechanics. Here, we 
show that it is in fact an advantage; distributions P{i) estimated from more about the system have a 
potential to describe the system better. The subjectivity becomes irrelevant in the thermodynamic 
limit owing to sharply peaked distribution P{C) dl])- In the current work, the sharply peaked limit 



arises naturally as limiting case of the maxEnt distribution (A — )> oo in Eq. 21 and C ~^ oo in 



Eq.|25|). 

We leave for future studies to investigate how constraint choice may dictate the predictive power 
of maxEnt based statistical mechanics. 
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V. SUPPLEMENTARY MATERIALS 
VI. NUMERICAL SIMULATIONS 

A harmonic spring consisting of two Lennard- Jones particles was immersed in a bath of 512 
Lennard-Jones particles in a cube of side 25A. NVT molecular dynamics simulations were run 
with NAMD (18). The CHARMM ()19p forcefield was used to describe the interaction between the 
Lennard-Jones particles and between the spring and the bath of particles. 

The spring constant for the harmonic oscillator was chosen to be A: = 0.5 kcal/mol-A^. The e 
parameter for the bath was set at —0.015 while the e parameter for the spring varied. We examined 
three different values of e = —5.5 and —10. The size parameter was set at r = 2.lA for the oscillator 
particles and r = 1.1 A for the bath particles. The systems were minimized for 2000 steps followed 
by an equilibration of Ins and a production run of 10ns. Configurations were stored every Ips. 
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VII. BEST FIT DISTRIBUTIONS 



In order to fit Eq. 21 and Eq. 25 to the experimental data, one needs to determine the free 
parameters from the simulation. In the traditional canonical ensemble, the inverse temperature 
/3 of the harmonic oscillator will be estimated from its average energy. Here, we show how to 
estimate the free parameters from the simulation. It is non-trivial to measure the average effective 
temperature /3 or the average system entropy {S{C)) in a computer simulation. Yet, operationally. 



(/?> 



1 



-^P{r\(3)P{(3)drdp 



{-)pP{P)dp. 



\p{r,P)drd^ 



(30) 



In other words, constraining /3 is equivalent to constraining 



Similarly, we can show that 



constraining S is equivalent to constraining logr. Thus, we estimate (r^), (^), and (logr) from 



the simulation and then fit Eq. 21 with (r^) and (logr), Eq. 25 (r^) and (^), and Eq. ?? from 
(r^>, (^>, and (logr). 



